Accurate particle position measurement from images 
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The moment method is an image analysis technique for sub-pixel estimation of particle positions. 
The total error in the calculated particle position includes effects of pixel locking and random noise 
in each pixel. Pixel locking, also known as peak locking, is an artifact where calculated particle 
positions are concentrated at certain locations relative to pixel edges. We report simulations to 
gain an understanding of the sources of error and their dependence on parameters the experimenter 
can control. We suggest an algorithm, and we find optimal parameters an experimenter can use to 
minimize total error and pixel locking. Simulating a dusty plasma experiment, we find that a sub- 
pixel accuracy of 0.017 pixel or better can be attained. These results are also useful for improving 
particle position measurement and particle tracking velocimetry (PTV) using video microscopy, in 
fields including colloids, biology, and fluid mechanics. 

PACS numbers: 52.27.Lw, 82.70.Dd, 07.05.Pj, 47.80.-v, 87.64.Rr 



INTRODUCTION 

Measurement of particle positions from images is im- 
portant in many fields, including dusty plasmas pHH, col- 
loids 0J3 , fluid mechanics [I| , biology @ , and computer 
vision [7j. Particle positions are generally estimated as 
the center of a bright spot of an image. Velocities can also 
be calculated from images; two common methods for this 
are Particle- Tracking- Velocimetry (PTV) and Particle- 
Image- Velocimetry (PIV) . 

To measure particle positions, an experimenter be- 
gins with a bit-map image. As an example, in Fig. [T] 
we present portions of single video frame from a dusty 
plasma experiment. Each bright spot represents an 8 /an 
diameter polymer microsphere illuminated by a 0.633 /im 
helium-neon laser sheet and imaged by a video camera 
with a Nikon 105 mm micro lens and a bandpass optical 
filter to eliminate unwanted light. The lens was focused 
to generate a sharp image. The experimental setup is 
similar to Fig. 2 of [lj. Figure HJa) and a magnified view 
Fig. Hfb) show portions of a video frame recorded by a 
cooled 14-bit digital camera (pcol600) with a 7.4 /im 
pixel width and a linear response. It was operated at 30 
frames per second with an exposure time of 30 msec. We 
should mention that experimental images of particles will 
differ, depending on many factors including the type of 
camera. To illustrate this point, we present in Fig. [TJc) 
an enlarged view of a bright spot in a frame recorded by 
an analog camera with a nonlinear response correspond- 
ing to gamma = 0.6. (Some cameras are nonlinear with 
an output intensity proportional to the input luminance 
to the power gamma). 

In the images in Fig. Q] particles fill several pixels. 
This spot size may be due, in part, to diffraction by the 
particle as well as camera properties such as diffraction 
by the camera aperture [8] and imperfect lens focusing. 
The spot size cannot be explained merely by geometrical 
optics, because the small particle size and magnification 



would result in an image smaller than one pixel on the 
camera detector. 

Images have random noise in each pixel. This can arise 
because of fluctuations in the camera's sensor and its elec- 
tronics. Noise in the experimental image of Fig. [TJa) is 
shown in Fig. [2] as a histogram of the pixel intensity. The 
most prominent feature is the noise peak, corresponding 
to a large number of pixels that are relatively dark. This 
noise peak has an average value that we term the "back- 
ground intensity," I(, g . The noise peak generally depends 
only on the camera and the sensor temperature. 

After recording a bit-map image, the experimenter will 
then use a computer algorithm to measure the particle 
position. There are several methods to do this, including 
the moment method 0, [3, 13" Hi) which we will study 
in this paper. Other methods include fitting a bright 
spot in the image to a Gaussian [12] or polynomial , 
and simpler methods such as choosing the centroid as 
the particle center @, [ll|. In the moment method, the 
calculated particle position is 
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where is the position and Ik is the intensity of a pixel 
k. The result of Eq. (JXJ) is sometimes called the "center of 
mass" 11] . When the particle fills more than one pixel, 



this calculation yields an estimate of the particle posi- 
tion with sub-pixel accuracy. Because of the efficiency 
and accuracy of the moment method, it is widely used 
when analyzing large quantities of data, as might be pro- 
duced for example when using a video camera. Fitting 
methods, which are more computationally expensive, are 
often used as well [I2| . The centroid method is similar to 
the moment method except that the intensity Ik of each 
pixel is replaced with a constant @, [H| ■ 

One application of particle position measurements is 
the calculation of particle velocities using PTV. A ve- 
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locity can be calculated by subtracting the positions of 
the same particle in two different frames and dividing by 
the time interval between frames. This method differs 
from PIV [l3j], where velocities are calculated at regular 
gridpoints, not for specific particles. 

Errors in the calculated particle position arise from 
multiple sources, including random noise in each pixel 
and also from the finite spatial resolution of the pixels on 
a camera sensor. When an image is recorded by sampling 
it with a finite number of pixels, some information about 
the intensity profile is lost, and this can cause a type of 
systematic error known as pixel locking or peak locking. 
The total error in the calculated position will be due to 
a combination of these effects, not just random noise or 
pixel locking by itself. 

In this paper, we seek to minimize the total error, and 
doing this will require that we understand the contribu- 
tion of pixel locking. Our goal is to aid the experimenter 
in making optimal choices, in both hardware and soft- 
ware, to minimize the total error. 



PIXEL LOCKING 

Pixel locking, also known as peak locking, is an arti- 
fact where calculated particle positions tend to be con- 
centrated at certain favored locations relative to pixel 
edges, such as the center or edges of a pixel. It is differ- 
ent from random errors, which do not result in favored 
positions for particles. To understand pixel locking, con- 
sider a particle whose image fills only a single pixel. In 
this case, the sum in Eq. (flj would have only a single 
term, and the position would be assigned to the exact 
center of that pixel. If the particle's image instead fills 
two pixels with equal intensity, the position will be as- 
signed to the midpoint of a pixel edge. The pixel center 
and midpoints of pixel edges are examples of favored po- 
sitions that are found to occur even when the particle's 
image fills several pixels 14j . 

The scientific literature for pixel locking includes many 
papers where PIV is used to measure velocities. In the 
early 1990s, the PIV method was tested to demonstrate 
their sub-pixel accuracy for particles flowing along with a 



fluid [15|, [16j . For specific applications of PIV, pixel lock- 
ing has been studied by other authors as well 17-2(J. in 



comparison to PIV, the literature for PTV includes fewer 
studies of pixel locking, e.g. [bj, |2l|. Because of this, 



some users of PTV, including until recently the authors 
of this paper, were unaware of pixel locking and the prob- 
lems it can cause. In addition to PTV, computer vision is 
another important area where pixel loc king is recognized 
as a problem in measuring positions [7|, [22} |23 [ . 

To detect pixel locking, we use sub-pixel maps as a 
diagnostic tool. A sub-pixel map shows all the calcu- 
lated particle positions relative to pixel edges, and it is 
drawn in a small box having the size of one pixel. To 



prepare a sub-pixel map, we begin with a graph of calcu- 
lated positions of N particles, as illustrated in Fig.[3ja), 
then plot the fraction parts of these positions in the small 
box, yielding the sub-pixel map in Fig.[3]Jb). In Fig.[3][c) 
we present an actual sub-pixel map calculated from a 
bit-map image by an analog camera in a dusty plasma 
experiment. The signature of pixel locking can be iden- 
tified in general by concentrations of calculated particle 
positions at favored positions. These favored positions 
can vary, depending on both hardware and software, but 
they commonly include the center or edges of a pixel, as 
in Fig. |3fc). Sub-pixel maps are therefore very useful for 
detecting pixel locking. Other authors have used similar 
graphs, where the calculated positions have been binned 
and plotted as a histogram 
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MOMENT METHOD 

The algorithm we optimize in this paper, the moment 
method, has two main steps. The first step is the se- 
lection of pixels that belong to each particle in the im- 
age. The second step is the calculation of position as an 
intensity-weighted moment of pixel positions. 

In the first step, the selection of pixels, the user begins 
by choosing a threshold I t h- The gray-scale image is re- 
placed by a black-and-white image, where pixels brighter 
than I t h become black, and all others become white. The 
choice of the threshold is important for several reasons 
Q, as we will discuss later. Next, the boundaries for in- 
dividual particle images are determined. There are sev- 
eral algorithms for selecting boundaries. We have exam- 
ined several codes that use the moment method, and we 
found that the only difference is the algorithm for se- 
lecting boundaries. We will consider three algorithms, 
which we distinguish by the corresponding codes we will 
test. All three of these codes are well tested, and they 
generate reliable results from experimental images. In 
one algorithm, the boundary is selected to be a poly- 
gon that encloses only contiguous pixels brighter than 
the threshold, Fig. |4ja). _This algorithm is used in the 



freely available ImageJ [l_l| code. The other two algo- 
rithms select a boundary that is a rectangle. In Code A, 
the boundary is the smallest rectangle that encloses all 
the contiguous pixels above the threshold 24], Fig. Sfb). 
In Code K, the boundary is the smallest rectangle that 
encloses a special curved contour (25|. This curved con- 
tour is produced by a 2D contour-plotting routine, and 
it is drawn not as line segments around pixel edges but 
rather as a curve passing through various pixels. Within 
a pixel, the pixel center is assigned the value of the orig- 
inal pixel intensity, but other points within a pixel are 
assumed to have other intensities, which are calculated 
by 2D interpolation using four surrounding pixel centers. 
Then, the contour-plotting routine draws a curve by join- 
ing all points, with sub-pixel spacing, where the assumed 
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intensity is equal to the threshold, as shown in Fig. |4}c) 
with a dash line. In both Codes A and K, but not Im- 
age J, the boundary can enclose some pixels that are less 
intense than the threshold. 

In the second step, which is the same in all three codes 
we test, the particle positions are calculated as the mo- 
ment, i.e., as the intensity- weighted position of pixels. 
The moment can be calculated 0, 0, Gil using Eq. (pj . 
However, we will find it better to use a generalized form 
of the calculated particle position, 



^calc 



Xfe(/fc — hase) 
_k 

J2(h — hase) 
k 



(2) 



where the baseline value hase will be explained later. 
Note that the calculated particle position depends on the 
selection of pixels that are included in the summation in 
Eq. or Eq. ©. 



METHOD 
Synthetic images 

To test methods of measuring particle positions, we 
calculate position errors as compared to true positions in 
synthetic images. For this purpose we cannot use actual 
experimental images because the true position is gener- 
ally not known. Synthetic images allow us to vary the 
intensity and the size of a bright spot to find how errors 
depend on these parameters. 

Units used in this paper are pixel units for all distances 
including X CQ i c , spot size and errors. Intensities, 
including signal and noise, are specified in intensity value 
units, i.e., a dimensionless integer ranging, for example, 
from to 2 14 — 1 for a 14-bit camera. 

We prepare synthetic images that resemble an experi- 
mental image like Fig. Ufa). The synthetic images have 
a size of 64 x 64 pixels, with one bright spot per image. 
These images have three major attributes that we com- 
pute: the spot's true position, the spatial profile of the 
signal, and the noise. 

First, the bright spot's true position is located near the 
image center, but displaced in the x and y directions by 
a fraction of a pixel. This is done using random numbers 
with a uniform distribution (between and 1) so that 
the true positions are random and uniformly-distributed 
relative to pixel edges. Using these random positions 
avoids any sampling bias. 

Second, like other authors [3,[26j], we model the signal's 
spatial profile as a Gaussian 
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characterized by the spot radius r spot and the peak in- 
tensity Ip e ak ■ (This Gaussian is intended to approximate 
the actual spatial profile, which depends on factors such 
as the particle size, the camera's gamma, and lens defo- 
cusing.) To imitate the collection of light onto a square 
pixel, we integrate this smooth profile over each pixel's 
area. This yields the value I S i g k of the signal in pixel k, 
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where k x and k y are the coordinates of pixel k. Equa- 
tion combined with Eq. © can be evaluated effi- 
ciently using the error function erf. (After this step, each 
bright spot has the same total signal intensity Yl I sig k , 
which was typically 37 707 corresponding to the brightest 
spot in the experimental image Fig. QJa) . In the experi- 
ment, not every bright spot has the same total signal in- 
tensity because some particles are levitated slightly above 
or below the brightest part of the horizontal laser sheet.) 

Third, we calculate a noise value Inoise k which is dif- 
ferent for each pixel k. To simulate the experiment, 
Inoise k is chosen as a random intensity from the noise 
distribution of our digital camera, Fig. [3J which is cen- 
tered at an average intensity hg = 384. Finally, the 
intensity Ik in each pixel is calculated as the sum of the 
intensities of the signal and noise or a saturation value 
I sat, whichever is smaller, 



Ik — Mm[(I sig k + Inoise k), I sat 



(5) 



(3) 



We use I sa t = 2 14 — 1 to simulate the saturation inten- 
sity of a real camera with 14-bit resolution. Finally, we 
round Ik to an integer because cameras produce integer 
values for the intensity of each pixel. The result of this 
calculation is a TIFF image like Fig. [SJb) or[5jc). 

Here we only consider bright spots that are circular, as 
in Eq. ([3]). Although we do not simulate them here, we 
note that non-circular bright spots can be analyzed using 
the moment method, and they do occur in some experi- 
ments. Elliptical particles arise when using analog video 
cameras with a limited horizontal resolution, or when 
particles move rapidly during the exposure time. The 
latter effect can be diminished by rastering a laser beam 
rather than dispersing it into a constant sheet. Defocus- 
ing a lens can result in non-circular spots, as in Sec. VII. 



Errors in calculated particle positions 

In this paper, we are mainly interested in errors in cal- 
culated particle positions. In addition to errors in par- 
ticle position, the experimenter may also be concerned 
with errors in velocities and other quantities computed 
from particle positions, as discussed in the Appendix. 
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To characterize the error in calculated particle posi- 
tions, we use two diagnostics. First, we calculate sub- 
pixel maps, as described in Sec. II. Examining these sub- 
pixel maps qualitatively will reveal pixel locking, which 
is one source of error. Second, we characterize the total 
error, including both random errors and pixel locking, as 
the root-mean-square (rms) difference of true and calcu- 
lated positions, i.e., the rms error 

" 1 N 1 5 

m— 1 

(6) 

where m and N are the index and total number, respec- 
tively, of bright spots. While we can calculate the total 
error using Eq. ([6]), we cannot separately calculate the 
contributions from random errors and pixel locking. 

To achieve good statistics, we prepared over 370 000 
synthetic images, each with one bright spot. We used 
N = 5000 when calculating the rms error, and N — 
100 000 when calculating sub-pixel maps. All of these 
images have different random true positions for their 
bright spots, and the noise in each pixel is different in 
all images. 



Parameters 

To find a procedure for calculating position with mini- 
mal total error, we will test three different codes, and we 
will vary parameters corresponding to software and hard- 
ware adjustments that an experimenter can make. We 
will now list these adjustments. The experimenter can 
choose to focus the camera lens sharply, or defocus it to 
make the bright spots in the image appear larger and fill 
more pixels. As a second parameter, the experimenter 
can adjust the image intensity by varying the camera 
aperture, exposure time or illumination brightness. After 
recording images with the camera, the experimenter will 
then use software. Here, we test three moment method 
codes, as explained in Sec. III. After choosing a code, the 
experimenter can usually adjust two parameters in that 
code: the threshold used in the first step, and the base- 
line (if any) that is subtracted in the second step, as in 
Eq. Q. 

Thus, we are motivated to analyze the impact of the 
following four parameters that the experimenter must 
choose: focus, intensity, threshold, and baseline. We do 
this by varying the values of r spot (keeping the total sig- 
nal intensity Isig k as constant, as explained later), 
Ipeak, hh, and I base , respectively. We will vary each of 
these four parameters in Sec. V. We will also compare 
results from the three different codes. The outcome of 
this analysis will be a practical procedure, presented in 
Sec. VI, that the experimenter can use to minimize errors 
in calculated positions. 



RESULTS 
Threshold 



The first parameter we vary is the threshold. The ex- 
perimenter will first choose a coarse range of threshold 
so that it is not so low that noise is wrongly identified as 
particles and not so high that fainter particles are over- 
looked. Then, within this coarse range, a fine adjustment 
can be made to reduce error. Here, we consider the fine 
adjustment. 

Our results in Fig. [6] show that the total error gener- 
ally increases with threshold, and it also depends on the 
choice of a code. We calculate the total error as the rms 
error, using N — 5000 images and Eq. ([6]). Recall that the 
total error includes both random and pixel-locking errors. 
The total error generally increases with the threshold be- 
cause raising the threshold can eliminate pixels that have 
useful signal. 

The total error exhibits not only a general increase 
with threshold, but also an oscillation. This is seen in 
Fig. [6l where there are several oscillations superimposed 
on the general trend. We cannot dismiss these oscilla- 
tions as mere statistical fluctuations because we achieved 
good statistics by using 5000 particle positions. To iden- 
tify the cause of these oscillations, we tested how the 
boundaries that are selected in the first step depend on 
the threshold. The result of this test is shown in Fig. [7] 
as a table of the boundaries selected by ImageJ. When 
the threshold is increased slightly so that the bound- 
ary shrinks by one pixel, there is a discrete jump in the 
calculated particle position. As the threshold increases, 
there is a sequence of jumps, as the boundary becomes 
smaller, one pixel at a time. These jumps, in aggregate 
for many particles, lead to oscillations in the rms error 
as the threshold is varied, which is the phenomenon we 
term the "boundary effect" . 

To identify the role of pixel locking in the total error, 
we examine sub-pixel maps in Fig. [8j which reveal the 
importance of the threshold. For ImageJ, we provide sub- 
pixel maps, Fig.[5][a) and^b), that correspond to the two 
thresholds that yielded the minimum and maximum rms 
errors, respectively, in Fig.|6] We note that the signature 
of pixel locking is weaker, i.e., the sub-pixel map is more 
uniform, for the case of the low threshold, Fig.[8ja), that 
yields the lowest total error. Conversely, the signature 
of pixel locking is stronger, i.e., the sub-pixel map has 
strongly non-uniform features, for the higher threshold, 
Fig. IHJb). In general, reducing the threshold will reduce 
pixel locking. Other codes exhibit the same trend, but 
with a different appearance for the sub-pixel maps, as in 
Fig. He) and^d). 
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Spot radius 

To simulate an experimenter's slight defocusing of a 
camera lens, we varied the spot radius r spot in Fig. [9J We 
used the Gaussian profile of Eq. ([3]), keeping the spot's 
total signal intensity (Yl Isig k summed over all pixels) 
constant. In this way we mimic an experiment where a 
particle scatters the same finite number of photons into 
a camera lens regardless of how the lens is focused. (We 
did not simulate the ring-shaped bright spot that can 
occur for extreme defocusing.) Defocusing can happen 
when an experimenter purposefully chooses to defocus 
the lens for example to avoid saturating pixels; in other 
cases, defocusing is not intentional but instead simply 
unavoidable because particles are at different depths, as 
for example in colloidal suspensions 0] and 3D dusty 
plasma suspensions [21 1. 

Defocusing a lens during the experiment can actually 
be desirable. By distributing the signal over a larger 
number of pixels, the impact of a single pixel in the cal- 
culation of the particle's position is less, so that pixel 
locking becomes weaker. On the other hand, defocusing 
can reduce the signal in each pixel, so that the signal-to- 
noise ratio (SNR) in each pixel becomes worse. In other 
words, there can be a trade-off: defocusing can improve 
pixel locking at the expense of making random errors 
worse. In our results below we investigate this effect. 

We should mention that when discussing defocusing, 
we always refer to the experimenter's adjustment to the 
hardware when recording an image. Unlike some other 
methods Q, here we do not blur an image in software 
after it has been recorded by the hardware. 

The result in Fig. [S] reveals three ranges of the spot 
radius, where the second range is the most desirable. In 
the first range, with small spot radii (r spot < 0.8), the 
total error diminishes with radius because the spot in- 
cludes a saturated pixel. Saturated pixels are undesirable 
because they introduce wrong information for intensity 
into Eq. @. In the second range, with slightly larger 
spot radii (0.8 < r spot < 2.0), the total error is smallest. 
In the third range, with large spot radii (r spot > 2.0), 
the total error generally increases with r spot because the 
trade-off results in the undesirable outcome of the wors- 
ened SNR in each pixel having a stronger effect than the 
improved pixel locking due to defocusing. The optimal 
spot radius is somewhere in the second range, which for 
our parameters is approximately 0.8 - 2.0. We should em- 
phasize, however, that this range will vary depending on 
the experiment due to different cameras (with different 
noise levels, sensitivities and saturation levels), particle 
size, illumination, and working distance between parti- 
cles and lens. If the camera had a higher noise level, the 
errors in this third range would be larger and the exper- 
imenter would be unable to use much defocusing. On 
the other hand, if the illumination were brighter, then 



the entire curve in Fig. [9J would shift toward larger spot 
radii and the experimenter would be able to use more 
defocusing. 

In Fig. [9J we also note an oscillation, superimposed on 
the general trend, for 0.8 < r spot < 2.0. We attribute this 
oscillation, which was observed previously in experiments 
by Kading and Melzer 21|, to a boundary effect similar 



to the one described above. 



Intensity 

To simulate adjusting the illumination brightness, the 
exposure time, or the camera aperture, we varied I pea k in 
Fig. [TO] As a result, the total signal intensity ^sig k is 
varied, while r spot is kept constant. We note that Image J 
yields the smallest total error. 

The trend that would be expected for random errors 
only is a downward slope as the intensity is increased, 
due to an improving SNR in each pixel. This trend is 
indeed observed Fig. [TUl but only for some of the data, 
as indicated by solid curves. The opposite trend is also 
observed in Fig. [TUl as indicated by dashed curves; since 
this trend is opposite to what is expected for random 
errors only, we attribute it to pixel locking. We term this 
particular effect of pixel locking the "pedestal effect." 



Baseline 

The pedestal effect is the result of a non-optimal choice 
of the baseline. To illustrate this effect, in Fig.[TT]we have 
sketched the cross section of a bright spot. The portion 
of this cross section that lies within the boundary, de- 
fined by the threshold, is shown shaded. This portion is 
divided in Fig. [11] into two parts, above and below the 
threshold. We term the part below the threshold the 
"pedestal," Fig. [TTJ The contribution of the pedestal to 
the moment in Eq. ^ can be large, or small, depending 
on whether Ibase is small or large, respectively. In the 
extreme case of a very large pedestal that dominates the 
calculation of the particle position, the calculated parti- 
cle position will often fall near a pixel edge or midpoint, 
as it does in the case of a centroid, thereby contributing 
to severe pixel locking. We term this tendency toward 
severe pixel locking the "pedestal effect." Below, we will 
determine the best choice of Ibase in order to reduce the 
pedestal effect and the pixel-locking errors that it intro- 
duces to the calculated particle positions. 

To test the effect of the baseline that is chosen, in 
Fig. [12] we present the total error, calculated as the rms 
error, for three different baseline values. From Fig. [12] 
we see that the total error is reduced by using a larger 
baseline value. The best choice is Ibase = hh, because 
this results in the smallest total error. It also minimizes 
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pixel locking; the downward slope in Fig. [12] indicates 
that random errors dominate. 

Thus, we conclude that in the second step, when using 
Eq. the baseline should be chosen to be the same as 
the threshold that was used in the first step. This can 
be done most simply by subtracting the same threshold 
for every pixel in the image. Alternatively, a different 
baseline level Ibase k for each pixel could be subtracted 
in Eq. ([2]) , to account for a different background level for 
each pixel. The latter method is useful because it allows 
the experimenter to eliminate optical reflections due to 
room lights, for example. The experimenter can calculate 
all the Ibase k baseline values for the pixels as follows. 
First, the experimenter will use the camera to record a 
"dark- field" image, with the illumination turned off so 
that particles are not visible. To improve the statistics, 
the experimenter can record a series of dark-field images 
and average them, pixel-by-pixel, to reduce the effect of 
random noise. This will yield an intensity Idark k for 
each pixel. Second, the baseline for each pixel will be 
calculated as 



',se k 



— Idark k + {hh ~ hg)- 



(7) 



Here, Ib g can be calculated as the average of Idark k f° r 
pixels in the image. 

With an optimal choice of both threshold and baseline, 
one can achieve a sub-pixel map that shows no evidence 
of pixel locking, as seen in Fig. riBT a). This map was pre- 
pared using Image J, with a baseline equal to the thresh- 
old. This choice of a baseline minimizes the total error, as 
we learned above. The reason that choosing Ibase = hh 
minimizes the total error is now clear: it greatly reduces 
pixel locking, so that mainly errors from random noise re- 
main. To further demonstrate the usefulness of choosing 
a baseline equal to the threshold, compare Fig. Ufa) to 
Fig. I13f a). The former figure, which was prepared sim- 
ilarly except with no baseline subtraction, reveals some 
pixel locking, while the latter does not. 

An experimenter, when attempting to choose optimal 
parameters, will be unable to calculate the rms error, 
as we have done in Fig. [121 f° r example. This is be- 
cause the true positions of particles are generally un- 
known. The experimenter can, however, calculate sub- 
pixel maps, such as Fig. [T3l because these require only 
calculated positions. Comparing Fig. [Tffia) and ITBTb') . 
which were both calculated with Ibase = Ithi but with a 
different hh, we see that the signature of pixel locking 
depends on the threshold. 

We now find our best result by varying the threshold, 
in Fig. 1141 to minimize the rms error. The threshold 
is the last parameter to choose, assuming that the ex- 
perimenter has already: (1) established the illumination 
level, (2) chosen a camera with a given noise level, (3) de- 
focused the camera lens to avoid saturating pixels, and 
(4) planned to use a baseline Ibase — Ith- Noting that 
the rms error in Fig. 1141 has several minima, we identify 



an optimal threshold by choosing the lowest minimum. 
This yields our best result, an rms error of 0.017. These 
same parameters also virtually eliminate the signature of 
pixel locking in Fig. I13f a). An experimenter can iden- 
tify an optimal Ith similarly, but without calculating the 
rms error, by examining sub-pixel maps for various val- 
ues of Ith, and among the maps with weak pixel- locking 
signatures, choosing the one with the lowest value of hh- 



PRACTICAL PROCEDURE 

We present here a practical procedure for using the 
moment method that minimizes the total error, includ- 
ing both random errors and pixel locking. This practical 
procedure includes first the use of hardware to record im- 
ages and then the use of software to analyze them. Our 
software uses the moment method with baseline subtrac- 
tion as we tested above; there are also other well-tested 
analysis methods that experimenters may wish to con- 
sider mm. 

For the hardware that produces the image, one will 
choose a camera and make adjustments to the intensity 
and lens focusing. Choosing a camera with low noise will 
not only reduce random errors; it will also allow the use 
of a lower threshold which can improve pixel locking. In 
using the camera, the optimal choices of intensity and 
lens defocusing must be considered together. The inten- 
sity can be varied, for example, by adjusting the cam- 
era aperture, exposure time, or illumination level. To 
achieve a high SNR in each pixel, we adjust the inten- 
sity upward as high as possible without saturating pixels. 
Another way to improve SNR is pixel binning, which also 
increase frame rate, but at the expense of spatial resolu- 
tion 



27]. If additional intensity is available but pixels are 



saturated, the experimenter can defocus the lens to avoid 
saturating the brightest pixels. Defocusing the lens helps 
reduce pixel locking, but it can increase random errors 
by reducing the SNR in each pixel; therefore, defocusing 
beyond a certain point actually worsens the total error. 
The optimal lens defocusing will depend on parameters 
such as intensity, camera noise level, and number of cam- 
era bits, which vary from one experiment to another. For 
the parameters we simulated (see Fig. we found that 
the optimal spot radius was in the range 0.8 - 2.0, mea- 
sured as the Gaussian half-width. For other parameters, 
we can offer this general guidance: the optimal lens de- 
focusing will be determined by the need to achieve an 
adequate SNR in each pixel. Noisier cameras or weaker 
illumination will require less defocusing, while low-noise 
cameras and brighter illumination will allow more defo- 
cusing. The lens should generally be defocused at least 
enough to avoid saturating pixels. 

For the image analysis software, there are usually three 
important choices. First, we prefer a code that has as its 
first step the selection of a boundary that includes only 
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contiguous pixels above a threshold. The freely avail- 
able ImageJ code selects such a boundary. Second, if 
the boundary is selected as described above in the first 
step, then in the second step, using Eq. ©, the base- 
line should be chosen equal to the threshold, in order to 
reduce pixel locking. This can be done either by sub- 
tracting the same baseline value from every pixel in a 
single step, or by using Eq. ([7]) with dark-field images if 
the experimenter wishes to remove the effect of optical 
reflections for example. Third, the threshold should be 
chosen in a two-part process. To start, the experimenter 
should count the number of particles that are identified, 
and then choose a coarse range as explained in Sec. V A. 
Next, within this coarse range, sub-pixel maps should be 
calculated for various thresholds. In order to reduce both 
random and pixel-locking errors, the user should choose 
the lowest threshold that has a weak signature of pixel 
locking. 

The moment method can achieve very low errors in 
particle position measurement when it is used optimally. 
For the case we simulated, an rms error as small as 0.017 
is achievable by making optimal choices in the software. 
Even smaller errors could be attained if the intensity were 
brighter or the camera had less noise. 

Readers who wish to perform tests similar to ours may 
use our codes and images [281 ]. 

EXPERIMENTAL DEMONSTRATION 

To demonstrate the practical procedure above, we used 
it in an experiment. The results presented above, based 
on synthetic images, indicate that both total errors and 
pixel locking will be reduced if we follow the practical 
procedure. Using experimental images, one can detect 
the signature of pixel locking using sub-pixel maps. We 
describe next the hardware and software components of 
our experimental test. 

For the hardware, the experiment was similar to the 
one for Fig. QJa), including using the same 14-bit cam- 
era, except that we improved the experimental method 
by slightly defocusing the lens. A cropped portion of the 
800 x 600 pixels image Fig. [TST a) and a magnified view 
Fig. Efb) show that a bright spot fills more pixels than 
in Fig. [IJb) where the lens was sharply focused. Due to 
defocusing, the spots are slightly noncircular. Addition- 
ally, we binned 2x2 pixels. As a result of these changes, 
the total intensity of a bright spot is typically 39 240, 
as compared to 21000 (with a maximum of 37 707) for 
Fig. Ufa), an d the noise peak is shifted to a lower inten- 
sity. A further possible improvement in the hardware is 
using a more powerful laser, and we plan to do that in 
future experiments. 

For the software, we used ImageJ to identify particles 
from 100 experimental images. We excluded any iden- 
tified particles that filled only one single pixel. First, 



we chose a coarse range for the threshold by counting 
the number of identified particles as a function of the 
threshold, Fig. [16] We looked for a nearly flat portion, 
which is from 325 to 925 here, and we chose that as the 
coarse range. Next, we calculated particle positions using 
Eq. ©, along with Eq. ([7]) to calculate Ib ase k using an 
average of 2000 dark-field images. We repeated these cal- 
culations of particle positions for various thresholds, each 
time preparing a sub-pixel map. Finally, we will examine 
these sub-pixel maps to choose the lowest threshold that 
has a weak signature of pixel locking. 

In Fig. 1171 we present the sub-pixel map that results 
from following our practical procedure in panel (a). Ex- 
amining this sub-pixel map, we see that it has no obvious 
signature of pixel locking when viewed in its entirety. To 
search for signatures, we zoom into the lower left corner, 
Fig J17f b)-(h). There, we can identify an artifact of pixel 
locking: a concentration of calculated positions on pixel 
edges. Our practical procedure requires choosing the low- 
est threshold with a weak signature of pixel locking. For 
our results in Fig. (TTJ thresholds in the range 325 - 425 
have no identifiable signature, leading us to choose 325. 

We conclude that the signature of pixel locking is 
vastly improved by using our practical procedure. This 
conclusion is based on a comparison of the sub-pixel maps 
in Fig. fTTT a) and Fig.(3][c). The latter was prepared for a 
similar experiment but a different camera, illumination, 
and analysis method. The signature of pixel locking is 
profound in Fig. Etc) , but it is virtually undetectable in 
Fig.&Mc). 
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APPENDIX: ERRORS IN OTHER QUANTITIES 

Errors in the calculated particle positions can intro- 
duce errors in other quantities that are calculated from 
the positions. In PTV, velocities are calculated as v = 
(x2 — xi)/At, as discussed in Sec. I. Pixel locking can af- 
fect the velocity calculation greatly in experiments. For 
example, if pixel locking is so severe that most calculated 
positions are located only at pixel centers, then almost all 
particle velocities calculated in PTV will be quantized as 
an integer number of pixel widths per frame. These errors 
in calculating velocities can propagate to other calcula- 
tions. Velocity distribution functions f(v) can be badly 
affected, with noticeable peaks [l4| that are signatures of 
pixel locking. However, we have found that wave spec- 
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tra and velocity correlation functions are not affected so 
badly. 

While it is beyond the scope of this paper to completely 
characterize the errors in v or f(v), we can discuss the 
contributions to the total error in v. For PTV, the rms er- 



ror, 5v = ((Sxi + 5x2 — 26x16x2) I Ai ) , has two con- 
tributions, [6x\ 



5x2 2 )/ At 2 arising from the errors in 



position, and (— 26x16x2) / At arising from correlations 
in the two errors. If the calculated position had random 
errors only, the correlation 6x16x2 would be zero and the 
rms error in v would be minimized when the rms error in 
x is minimized. However, pixel-locking errors can have 
correlations, which will vary depending on the velocities, 
and these will affect 5v in a way that is difficult to pre- 
dict. 

Aside from these quantities, which are calculated from 
velocities, experimenters often calculate other quantities 
from the position itself. The mean-square displacement 
(MSD), which is used to measure diffusion, is calculated 
from position. Particle position errors can cause the MSD 
to be exaggerated significantly at small times when the 
displacement is small, but not at large times when the 
displacement is large |f|. Another use of particle posi- 
tions is the study of structure 29|, 3(| . While we have not 
analyzed the sensitivity of structural analysis methods to 
particle position errors, we expect that calculations that 
are sensitive to small changes in interparticle distances, 
such as Voronoi maps for detecting defects, will be more 
affected than correlation functions that use data over a 
wide range of distances. 
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FIG. 1: Experimental bit-map images of a monolayer sus- 
pension of microspheres in a dusty plasma. Each bright spot 
corresponds to one particle. Here, (a) is 1/12 of the origi- 
nal image from a digital camera and (b) is a magnified view, 
showing that a bright spot fills several pixels, while in (c) from 
an analog camera a bright spot fills about 5x5 pixels. Spot 
size depends on such factors as camera type and focusing. A 
particle's position is calculated as the bright spot's center; 
errors in this calculation are the topic of this paper. 



FIG. 2: Histogram of intensity values of pixels in the original 
experimental image of Fig. Ola). The inset shows the same 
data with a logarithmic scale. The prominent peak, centered 
at It g , is due to noise in the camera. 



FIG. 3: Illustration of the method for calculating a sub-pixel 
map. First, a 10 x 10 pixel bit-map image (not shown here) is 
analyzed to yield a map (a) of particle positions. Second, the 
same positions are plotted relative to pixel edges in (b); these 
values are the fraction parts of the calculated positions, (c) 
An example sub-pixel map of N = 617 particles, calculated 
from an experimental image (full view of Fig. [He)), reveals 
pixel locking as a tendency of calculated positions to be con- 
centrated at favored positions including the center and edges 
of pixels. 



FIG. 4: Illustration of boundaries. In algorithms for calcu- 
lating particle positions from a bit-map image, the first step 
is selecting the contiguous pixels to be used, as defined by a 
boundary (solid white line) that encloses them. The codes 
tested here differ only in the way they select boundaries, (a) 
In ImageJ, only contiguous pixels above a threshold are in- 
cluded in the boundary. Code A (b) and Code K (c) use 
boundaries that are the smallest rectangles that enclose: all 
the contiguous pixels above the threshold in Code A, or the 
dashed contour produced by a 2D contour-plotting routine in 
Code K. 



FIG. 5: Magnified images of bright spots, (a) Experimental 
image from a digital video camera. (b),(c) Synthetic images, 
with a Gaussian profile centered on a known true position, 
here with two different spot radii. In generating synthetic 
images, we first choose the true position randomly, and then 
calculate the intensity of each pixel using Eq. (J3J) so that it 
includes both signal and noise. 



FIG. 6: The rms error of calculated positions as a function of 
the threshold I t h- In general, errors increase with threshold, 
and superimposed on this increase is an oscillation. The rms 
errors are always calculated as in Eq. ^ using N = 5000. 
(Here, r 3po t = 1.5 pixel units, I pea k = 5334 intensity value 
units, corresponding to a total signal intensity ^ 7 S j s k = 
37 707. Also, I base = 0.) 



FIG. 7: Cause of oscillations. Boundaries, selected in the 
first step of Image J, enclose fewer pixels as the threshold is 
increased. Removing one pixel from the boundary causes a 
discrete jump in the calculated particle position in Eq. (T2|. 
As the threshold increases, there is a sequence of jumps, as 
the boundary becomes smaller, one pixel at a time. These 
jumps, in aggregate for many particles, lead to oscillations in 
the rms error as the threshold is varied, a phenomenon we 
term the boundary effect. The three columns correspond to 
three different true positions. 



FIG. 8: Sub-pixel maps for N = 100 000 randomly dis- 
tributed true positions. The signature of pixel locking is gen- 
erally more severe for higher thresholds. (Here, r apot = 1.5, 
Ipeak = 5334, and I base = 0.) 



FIG. 9: Simulation of slight lens defocusing. The optimal 
range of spot size lies between two other ranges: for very small 
r 3p ot, errors worsen due to pixel saturation; for very large 
r spo t, they worsen due to random errors. For our parameters, 
these two ranges are for r apot < 0.8 and r spot > 2.0, respec- 
tively. Oscillations in the optimal range arise from a boundary 
effect. (Here, 7 th = 1000, I base = 0, and £)iai fl k = 37 707.) 



FIG. 10: The rms error as the intensity is varied, to simulate 
adjusting the illumination brightness, the exposure time or 
the camera aperture. The main trend is that the error de- 
creases with increasing intensity due to an improved signal- 
to-noise ratio (SNR), as indicated by solid curves; the op- 
posite trend, indicated by dashed curves, is attributed to a 
pixel-locking effect that we term the pedestal effect. (Here, 
r sp ot = 1.5, I t h = 740, and I base = 0.) 



FIG. 11: Cross section of a bright spot, illustrating the 
"pedestal." Pixels brighter than the threshold identify the 
boundary for ImageJ in the first step. In the second step, 
both shaded portions contribute to the calculated particle po- 
sition if hase = 0, i.e., if no baseline is subtracted in Eq. 
The lower shaded portion, marked "pedestal," can heavily in- 
fluence the calculated particle position. The pedestal can be 
reduced by choosing hase ~ hg, or eliminated altogether by 
choosing I base = I th . 



FIG. 12: Test of different baselines. The best choice to mini- 
mize rms error is subtracting a baseline equal to the threshold 
Ith in Eq. ([2}. (We used ImageJ, and r spo t = 1-5, Ith = 740, 
and I bg = 384.) 



FIG. 13: Sub-pixel maps, using a baseline hase = Ith for two 
different thresholds (a) I th = 1150 and (b) I th = 2950. Com- 
paring these panels shows that the signature of pixel locking 
can be virtually eliminated, as in (a), by making the best 
choice of threshold as well as choosing I base = I t h- (Here, we 
used the same 100 000 images as in Fig. [8]) 



11 



FIG. 14: Total error, using a baseline Ibaae = Ith- Comparing 
to Fig. 6 where Ibaae = 0, errors have been reduced. The low- 
est rms error that can be achieved with these images is 0.017, 
using the same optimal choice of parameters as in Fig. I13f a~). 
(We used the same 5000 images as in Fig. [6] Here and in 
Fig- El we used ImageJ.) 



FIG. 15: Experimental bit-map images of a monolayer sus- 
pension of microspheres in a dusty plasma. Here, (a) is 1/12 
of the original image and (b) is a magnified view. A bright 
spot fills about 5x5 pixels. Compared to Fig.QJa), the hard- 
ware was improved by slight lens defocusing. 



FIG. 16: Choosing the coarse range of threshold using ex- 
perimental images. Counting the particles identified in 100 
images, we choose the nearly flat portion 325 < I t h < 925 as 
the coarse range. Outside this coarse range, many false parti- 
cles appear at lower I t h due to noise, while many true particles 
are missed at higher I t h- Labels a-h identify thresholds used 

in Fig. era 



FIG. 17: Experimental sub-pixel maps for different thresholds 
within the coarse range. Here, (a) is an entire map, and 
(b)-(h) show the lower left corner. We choose the lowest I t h 
with a weak signature of pixel locking, 325. The signature 
is stronger for iy, > 525, with a concentration of calculated 
positions on pixel edges. Vastly better than Fig. EJc), there 
is no obvious signature of pixel locking for I t h < 525. (Here, 
we used ImageJ with Ib ase k calculated from Eq. (O and a 
dark-field image.) 
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